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ABSTRACT 


The  objective  of  this  research  program  is  to  improve  the  capability  to  predict  the  seismic  source 
characteristics  of  underground  explosions  in  rock.  This  is  being  accomplished  by  development  of  improved 
dynamic  failure  models  that  are  constrained  by  a  large  unique  data  set  of  near-field  waveforms  and 
parametric  data  from  both  U.S.  tests  and  historic  Soviet  explosions  at  the  Degelen  Test  Site.  In  addition,  we 
are  analyzing  regional  seismic  data  along  a  seismic  line  located  north  of  the  Degelen  Test  Site  that  recorded 
data  at  9  stations  spaced  approximately  evenly  from  the  test  site  to  a  distance  of  about  100  km.  This  project 
is  a  collaborative  effort  between  SAIC  and  the  Russian  Institute  for  the  Dynamics  of  the  Geospheres  (IDG). 

IDG  is  in  the  process  of  digitizing  data  from  25  nuclear  tests  at  Degelen.  The  complete  data  set  consists  of 
81  near-field  waveforms  recorded  underground  at  shot  depth  and  124  near-regional  seismic  records.  Most 
seismic  records  include  both  a  vertical  and  radial  waveform.  This  data  set  provides  a  rare  opportunity  to 
observe  and  model  the  seismic  wavefield  of  the  explosions  as  they  evolve  from  the  near  field  of  the 
explosion  out  to  regional  distances.  To  date,  IDG  has  digitized  27  near-field  waveforms  and  122  seismic 
waveforms.  We  are  in  the  process  of  defining  the  optimum  procedures  to  model  this  data.  The  goal  is  to 
develop  material  models  that  are  consistent  with  the  data  and  have  a  realistic  physical  basis.  Work  to  date  at 
SAIC  has  focused  on  implementation  and  testing  of  improved  numerical  modeling  procedures  and 
simulation  of  near-regional  data.  We  are  testing  acoustic  fluidization  as  a  physical  mechanism  for  strength 
reduction  in  nonlinear  explosion  simulations.  Near  regional  data  is  being  modeled  using  wave  number 
integration.  We  are  also  comparing  the  Degelen  data  with  similar  data  from  NTS  explosions.  The  near-field 
Degelen  data  received  to  date  is  lower  in  amplitude  than  would  be  predicted  based  on  NTS  experience. 
Although  both  the  NTS  explosions  being  reviewed  and  the  Degelen  explosions  were  in  granite,  the  material 
properties  of  the  rock  are  significantly  different.  The  most  significant  difference  is  higher  porosity  in  the 
Degelen  granite.  The  data  set  of  near-field  and  regional  waveforms  will  be  delivered  to  the  Center  for 
Monitoring  Research  upon  completion  of  this  project. 

KEY  WORDS:  rock  mechanics,  Degelen,  explosion,  numerical  modeling 
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OBJECTIVE 


The  main  objective  of  this  research  program  is  to  improve  the  capability  to  predict  the  seismic  source 
characteristics  of  underground  explosions  in  rock.  This  is  being  accomplished  by  development  of  improved 
dynamic  failure  models  for  hardrock  that  are  constrained  by  a  large  unique  data  set  of  near-field  waveforms 
and  parametric  data  from  both  U.S.  tests  and  historic  Soviet  explosions  at  the  Degelen  Test  Site.  In 
addition,  we  are  analyzing  regional  seismic  data  along  a  seismic  line  located  north  of  the  Degelen  Test  Site 
that  recorded  data  at  9  stations  spaced  approximately  evenly  from  the  test  site  to  a  distance  of  about  100 
km.  This  project  is  a  collaborative  effort  between  SAIC  and  the  Russian  Institute  for  the  Dynamics  of  the 
Geospheres  (IDG). 

RESEARCH  ACCOMPLISHED 


Introduction 


Empirical  and  numerical  models  of  explosion  sources  do  a  fairly  good  job  of  matching  observed  seismic 
signals,  but  the  physical  basis  for  the  explosion  source  is  still  not  well  understood.  In  particular,  numerical 
models  of  explosion  sources  developed  using  laboratory  measurements  of  rock  properties  predict  particle 
velocity  pulse  widths  that  are  roughly  a  factor  of  3-5  smaller  than  those  observed  near  field  in  hardrock 
events  such  as  PILEDRIVER.  The  basic  problem  is  that  the  strength  of  the  rock  measured  in  the  laboratory 
is  much  larger  than  the  apparent  strength  of  the  rock  as  determined  from  the  near  field  ground  motion. 

A  number  of  solutions  to  these  problems  have  been  proposed  over  the  years,  including  the  effective  stress 
model  (Rimer,  et  al,  1984),  and  various  types  of  damage  models.  These  models  all  have  the  characteristic 
that  the  material  strength  is  reduced  dynamically  to  a  very  low  level  when  the  rock  fails.  The  effective 
stress  model  says  that  the  weakness  comes  from  water  within  the  rock  matrix,  and  that  the  broken  rock  in 
effect  floats  on  water  that  is  squeezed  out  of  pores  or  fractures  when  the  rock  fails.  Although  there  are 
questions  about  the  realism  of  this  physical  model,  it  works  fairly  well  to  explain  the  near  field  waveforms. 
Figure  1  shows  a  comparison  of  the  PILEDRIVER  waveforms  with  waveforms  calculated  using  the 
effective  stress  model.  The  agreement  is  quite  good,  particularly  at  the  closer  two  stations.  Furthermore, 
when  the  PILEDRIVER  solution  was  scaled  to  the  appropriate  yield  and  compared  with  other  US 
explosions  in  granite  (Hardhat  and  Shoal),  agreement  with  the  observed  waveforms  was  also  quite  good 
(Stevens,  et  al,  1986). 


Figure  1.  Particle  velocity  measurements  at  working  point  depth  from  PILEDRIVER  compared  to 
numerical  simulation  pile570  using  the  effective  stress  law. 
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Under  a  previous  DTRA  contract,  SAIC  (at  the  time  Maxwell  Technologies)  worked  together  with  the  IDG 
and  the  University  of  Southern  California  (USC)  to  develop  improved  material  models.  IDG  provided 
extensive  measurements  of  material  properties  close  to  nuclear  and  chemical  explosions  both  before  and 
after  the  explosions  were  detonated  (Rimer,  et  al,  1998).  In  addition,  we  implemented  the  micro¬ 
mechanical  damage  mechanics  model  for  the  growth  and  coalescence  of  cracks  in  brittle  rocks  under 
compressive  loading  which  was  developed  by  Prof.  Charles  Sammis  of  USC  into  SAIC’s  nonlinear  finite 
difference  codes  and  then  used  this  model  to  simulate  the  observed  explosion  damage  and  a  small  set  of 
near  field  waveforms  that  were  also  provided  by  IDG.  The  results  of  this  work  are  discussed  in  detail  in  the 
final  report  (Rimer  et  al,  1999).  Since  the  Sammis  damage  model  does  not  predict  what  happens  to  the  rock 
after  unstable  brittle  failure  occurs,  the  calculations  introduced  a  friction  law  model  post-failure.  However, 
realistic  values  of  friction  did  not  provide  enough  strength  reduction  to  match  the  near  field  data.  We  were 
more  successful  in  matching  the  data  by  dropping  the  coefficient  of  friction  to  very  low  values  (as  low  as 
0.02)  for  rubbleized  rock,  but  this  again  leaves  the  question  of  what  physical  mechanism  could  be 
responsible  for  these  very  low  values  and  corresponding  low  strength. 

A  possible  answer  initially  proposed  by  Melosh  (1979)  is  “acoustic  fluidization’’.  The  physics  behind  this 
mechanism  is  that  there  is  a  complex  dynamic  acoustic  wavefield  that  causes  high  frequency  vibrations  in 
the  broken  rock.  These  vibrations  cause  rapidly  changing  regions  of  high  and  low  normal  stress,  and 
remove  the  frictional  normal  stress  from  parts  of  the  rock  as  it  moves.  Consequently  parts  of  the  rock  are 
not  confined  by  the  frictional  stress  and  in  effect  have  much  lower  strength.  Acoustic  fluidization  has  been 
used  to  explain  other  phenomena  such  as  landslides  and  craters  (Melosh  and  Ivanov,  1999),  which  have 
been  similarly  difficult  to  explain  because  of  anomalously  low  apparent  friction.  Initial  efforts  to  include 
acoustic  fluidization  in  our  numerical  models  are  described  later  in  this  paper. 

Degelen  Near  Field  and  Near  Regional  Seismic  Data 

The  numerical  modeling  component  of  this  project  is  being  constrained  by  a  much  better  data  set  than  has 
been  available  in  the  past.  Near  field  waveforms  are  only  available  from  a  small  number  of  U.S.  nuclear 
tests,  and  until  recently,  none  have  been  available  from  the  testing  program  of  the  former  Soviet  Union. 
IDG  has  near  field  records  from  a  number  of  nuclear  explosions  at  the  Degelen  test  site  that  are  being 
digitized  for  this  project.  IDG  will  also  be  providing  near  source  material  properties  measurements  from 
these  events.  In  addition,  IDG  has  data  from  a  seismic  line  located  north  of  the  Degelen  test  site  that  was 
maintained  with  consistent  instrumentation  for  many  years  during  the  Soviet  testing  program.  Figure  2 
shows  a  map  illustrating  the  location  of  the  test  sites  and  the  seismic  stations. 


Figure  2.  Map  showing  the  locations  of  the  former  Soviet  Degelen  and  Balapan  test  sites,  faults,  and 
seismic  stations  in  the  region. 
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The  9  seismic  stations  are  spread  out  at  approximately  even  intervals  from  the  Degelen  test  site  out  to  a 
distance  of  about  100  km.  IDG  is  also  digitizing  seismic  data  from  these  stations  for  all  of  the  events  that 
have  near  field  records.  This  provides  a  rare  opportunity  to  observe  and  model  the  seismic  wavefield  of  the 
explosions  as  they  evolve  from  the  near  field  of  the  explosion  out  to  regional  distances.  The  data  that  IDG 
has  identified  for  digitization,  and  the  data  digitized  as  of  this  writing  are  listed  in  Table  1.  Figure  3  shows 
waveforms  from  one  of  these  explosions. 


Table  1.  Degelen  events  with  near  field  and/or  seismic  records,  and  waveforms  digitized  to  date.  (The  total 
number  of  digitized  records  including  multiple  recordings  at  some  stations,  are  shown  in  parentheses. 
Explosion  yield  and  ISC  mb  are  shown  for  events  digitized  to  date.) 


Date 

Yield 

mb 

Number  of  Waveforms 

Distance 

Number  of 

Waveforms 

Distance 

(Kt) 

near  field  records  Digitized 

range,  m 

seismic  records 

Digitized 

range,  km 

1988/04/22 

2.3 

4.9 

2 

2 

69-183 

3 

3 

57-81 

1988/11/23 

19 

5.4 

4 

3 

280-475 

9 

9 

14-83 

1987/07/17 

78 

5.8 

14 

14 

170-900 

7 

7 

14.5-84 

1987/10/16 

1.1 

4.6 

4 

4 

55-132 

8 

8  (27) 

20-76 

1989/10/04 

1.8 

4.7 

4 

4(5) 

45-285 

8 

8(24) 

16-85 

1981/07/17 

9.3 

5.1 

3 

3 

115-310 

8 

8  (27) 

15-80 

1965/02/04 

.001-20 

- 

4 

4 

147-750 

9 

peaks 

14-83 

1964/05/16 

20-150 

6.2 

4 

4 

150-600 

9 

peaks 

13-80 

1966/03/20 

4 

300-600 

9 

15-84 

1966/02/13 

4 

350-600 

9 

14-80 

1973/12/31 

4 

110-230 

1974/12/16 

4 

100-200 

1978/03/26 

5 

76-270 

1980/06/25 

3 

155-310 

1982/12/25 

4 

90-190 

1984/10/18 

4 

35-110 

1968/07/12 

2 

80-90 

1968/12/18 

3 

200-510 

1970/06/28 

2 

200-240 

1980/09/25 

3 

100-160 

1971/12/15 

9 

7.5-85 

1987/05/06 

9 

13-83 

1987/12/20 

9 

13-83 

1988/10/18 

9 

11-77 

1989/02/17 

9 

10-77 
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Figure  3.  Near-field  (left)  and  near  regional  seismic  (right)  waveforms  from  the  1987/07/17  event. 

The  near  field  records  in  Figure  3  show  the  evolution  of  the  waveform  from  the  nonlinear  to  linear  regions. 
Unfortunately,  the  absolute  times  are  not  known  for  the  records.  The  near  regional  records  show  the 
evolution  of  the  waveform  from  14  to  83  km.  A  strong  Rg  phase  is  present  in  several  seismograms  that 
persists  to  quite  a  long  distance.  Some  of  the  records  end  before  the  start  of  the  Rg  phase  and  therefore  do 
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not  show  it.  Figure  4  shows  a  comparison  of  synthetic  and  observed  seismograms  for  three  of  the  near 
regional  waveforms  from  the  1987/07/17  event.  The  synthetic  seismograms  were  constructed  using 
wavenumber  integration  (Luco  and  Apsel,  1983)  using  the  East  Kazakh  structure  from  Stevens  (1986).  The 
synthetic  seismograms  were  low-pass  filtered  at  2  Hz.  The  persistence  of  Rg  calls  into  question  the 
explanation  of  Rg  scattering  as  the  source  of  Lg  since  that  mechanism  requires  most  of  Rg  to  scatter  into 
Lg  within  a  few  kilometers  of  the  source. 

Figure  5  shows  a  comparison  between  peak  velocity  measurements  from  the  first  three  explosions  listed  in 
Table  1  and  peak  velocity  measurements  from  other  explosions  in  granite.  The  Degelen  velocities  are  near 
the  lower  bound  of  the  other  velocity  measurements. 
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Figure  4.  Synthetic  and  observed  seismograms  for  the  1987/07/17  event.  Synthetics  were  created  using 
wavenumber  integration  and  were  low  pass  filtered  at  2  Hz. 


M5  Jua 

u 


14  in 


77,  El 

1  a  ki 

2  3  KI 


i  fPliSh  CATA 
.  ralE  4lE 
H  I4HL 
HiPrrr.vTi 


5 


■;  , 


Tv. 


ir 

fe  ■  h  1 1 .  i 


Figure  5.  Peak  velocity  vs.  scaled  range  for  the  first  three  Degelen  explosions  (left)  and  for  historic  data 
from  U.S.  and  French  explosions.  The  solid  line  is  the  prediction  from  the  PIFEDRIVER  simulation 
discussed  in  the  next  section  and  shown  earlier  in  Figure  1.  The  Degelen  velocities  are  near  the  lower 
bound  of  the  historic  data  set. 
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Numerical  Modeling  of  Explosions 

The  parameters  used  in  numerical  simulations  of  underground  nuclear  explosions  are  constrained  by 
laboratory  material  properties  tests  and  by  direct  observations  of  ground  motion  from  underground 
explosions.  Laboratory  measurements  of  strength  for  brittle  hardrocks  seem  to  be  very  inconsistent  with 
insitu  strength  as  inferred  by  modeling  of  underground  explosions.  In  particular,  finite  difference 
calculations  of  ground  motion  in  granite,  made  using  the  laboratory  measurements  of  shear  strength,  have 
invariably  given  much  narrower  particle  velocity  pulses  and  much  smaller  displacements  than  those 
measured  in  the  field.  Constitutive  models  have  been  developed  (Rimer  and  Lie,  1982,  Rimer,  et  al.,  1984) 
which  attribute  this  weaker  behavior  of  insitu  granite  and  other  rocks  under  explosive  loading  to  ground- 
motion-induced  rock  damage  or  pore  fluid  pressure  increases.  Comparisons  have  already  been  made  in 
Figure  1  between  the  particle  velocity  measurements  at  working  point  depth  from  the  62  kt  PILEDRIVER 
event  and  the  results  of  the  pile570  numerical  simulation,  made  using  the  effective  stress  model  discussed 
in  Rimer,  et  al  (1984).  Numerical  results  from  pile570  are  also  in  good  agreement  with  the  measured  cavity 
radius  and  with  the  estimated  seismic  source  function  given  in  Table  2. 


Table  2.  United  States  Explosions  in  Granite 


Explosion 

Yield  (kT) 

Measured  Cavity  Radius  (m) 

Calculated  (570)  Cavity  Radius  (m) 

PILEDRIVER 

62.0 

40.1/44.5 

42.5 

HARDHAT 

5.9 

19.4 

19.4 

SHOAL 

12.5 

26.8 

24.9 

For  this  simulation,  effective  stress  model  parameters  were  calibrated  to  best  match  the  velocity  peak  and 
pulse  shape  at  the  closest-in  gauge  station  (B-SL).  This  required  a  rapid  buildup  of  pore  pressure  during 
the  loading,  leading  to  a  large  strength  reduction  very  near  the  propagating  shock  front.  Note  that  the  pulse 
shapes  at  all  four  stations  are  rather  consistent,  with  all  including  a  long  duration  shallow  negative  velocity 
pulse.  The  peak  velocities  at  the  two  smallest  ranges  agree  very  well  with  the  simulation,  but  the  measured 
peaks  at  the  two  larger  ranges  are  a  factor  of  two  or  more  lower  than  the  calculated  peaks.  The  two  closest 
PILEDRIVER  gauges  were  located  along  roughly  a  180  degree  different  azimuth  than  the  other  more 
distant  gauges.  However,  a  connection  between  possible  site  anisotropies  and  the  measured  ground  motions 
has  never  been  established.  It  should  be  noted  that  the  5.9  kt  HARDHAT  event,  which  was  detonated  in 
similar  rock,  near  the  later  PILEDRIVER  event,  but  at  a  shallower  depth  of  burial,  gave  particle  velocity 
pulses  that  were  more  similar  to  the  PILEDRIVER  pulses  at  the  larger  ranges.  (See  Rimer,  Stevens,  and 
Day,  1987,  for  a  comparison  between  the  results  of  the  pile570  simulation  and  the  HARDHAT 
measurements.) 

The  constitutive  models  used  in  these  simulations  are  phenomenological  in  nature  and  do  not  explicitly 
account  for  the  dynamic  response  of  the  insitu  fractures  in  the  crystalline  rock.  The  Sammis  micro¬ 
mechanical  damage  model  for  brittle  rocks  under  compressive  loading  (see  Ashby  and  Sammis,  1990)  has 
been  incorporated  into  our  numerical  simulation  codes.  This  model  introduces  a  damage  parameter,  related 
to  the  increase  in  flaw  size  from  its  pre-shot  average  value.  Damage  accumulates  as  the  flaws  extend 
during  the  compressive  loading  and  reaches  some  maximum  value  at  which  the  rock  fails  unstably.  Since 
unstable  compressive  failure  of  a  rock  element  is  calculated  using  this  model  to  occur  relatively  early  in  the 
dynamic  motions  of  interest  here,  i.e.,  usually  near  the  propagating  shock  front,  additional  modeling  was 
incorporated  to  complete  the  description  of  the  stress  field  after  this  failure.  Limiting  the  magnitude  of 
deviatoric  stresses  in  a  failed  rock  element  through  the  use  of  a  standard  friction  law  was  shown  to  not 
provide  sufficient  strength  reduction  to  simulate  the  ground  motion  measurements.  In  particular,  calculated 
particle  velocity  time  histories  were  still  much  narrower  than  the  measured  pulses  for  all  reasonable  choices 
of  initial  flaw  size. 

The  additional  strength  reduction  required  to  sufficiently  broaden  the  particle  velocity  pulses  was  obtained 
by  using  a  shear  damage  model  originally  developed  for  soft  rocks  such  as  tuff,  described  in  Rimer  and 
Proffer  (1991).  As  discussed  in  Rimer  et  al  (1999),  this  shear  damage  model  was  applied  here  only  for  rock 
elements  that  had  experienced  the  onset  of  unstable  compressive  failure.  Thus,  a  rock  element  was  allowed 
to  undergo  significant  damage  before  application  of  this  treatment.  The  post-failure  shear  damage  model 
performs  an  interpolation  between  the  standard  coefficient  of  friction  of  0.60  and  much  lower  “effective 
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friction”  values  of  0.02-0.20  in  the  friction  law  used  to  limit  deviatoric  stresses.  This  linear  interpolation  is 
based  on  the  maximum  shear  strain  experienced  by  the  failed  rock  element. 

A  series  of  calculations  were  made  using  the  Sammis  modeling,  quantifying  the  effect  of  model  parameters 
on  particle  velocity  pulse  shapes,  cavity  radius,  and  RDP.  Calculated  particle  velocity  pulse  widths  were 
sensitive  to  the  choice  of  shear  strain  magnitude  required  for  full  strength  reduction,  with  lower  values  of 
this  parameter  resulting  in  longer  pulse  duration  and  larger  RDP.  The  simulation  which  best  fits  all  of  the 
PILEDRIVER  ground  motion  data  is  Run  PD  10.  Comparisons  with  measured  particle  velocities  are  shown 
in  Figure  6.  In  contrast  with  the  results  of  pile570.  Figure  1,  made  with  the  effective  stress  model,  PD10 
provides  a  much  better  fit  to  the  PILEDRIVER  data  at  the  two  larger  ranges,  while  underestimating  peak 
velocities  at  the  closest  ranges.  Subsequent  analyses  showed  that  the  timing  of  the  strength  reductions  for 
the  two  models  were  somewhat  different,  with  the  effective  stress  model  providing  an  earlier  reduction  than 
the  present  damage  model.  However,  calculated  cavity  radius  and  RDP  with  the  new  model  are  in  as  good 
agreement  with  the  measurements  as  were  those  with  the  effective  stress  model.  It  is  important  to 
emphasize  that  while  the  shear  damage  model  (or  some  other  post-failure  strength  reduction  model)  is 
crucial  to  successful  simulation  of  the  ground  motion  data,  it  is  the  micro-mechanical  damage  mechanics 
model  that  primarily  determined  the  size  of  the  central  core  of  low  strength  rock  around  the  cavity. 


Time  (seconds) 


Time  (seconds) 


Figure  6.  Calculated  particle  velocity  pulses  at  four  ranges  for  Run  PD10  and  PILEDRIVER 
measurements 

Acoustic  Fluidization 


A  possible  physical  mechanism  for  the  strength  reduction  and  very  low  coefficients  of  friction  discussed  in 
the  previous  section  is  the  “acoustic  fluidization”  concept  proposed  by  Melosh  (1979)  to  explain  the  low 
strengths  (or  very  low  angles  of  internal  friction)  apparent  in  a  number  of  geologic  processes,  such  as 
seismic  faulting,  impact  crater  slumping,  and  long  runout  landslides.  The  main  concept  of  this  proposed 
mechanism  is  that  “sufficiently  strong  acoustic  waves  in  rock  debris  can  momentarily  relieve  the  static 
overburden  in  limited  regions  of  the  rock  mass,  thus  allowing  sliding  to  occur  in  the  unloaded  region.  If 
this  happens  frequently  enough,  flow  of  the  entire  rock  mass  results.”  In  terms  of  the  explosively-induced 
ground  motions  of  interest  here,  only  the  mass  of  failed  or  fractured  rock  would  have  the  potential  for 
acoustic  fluidization.  The  dynamic  fracturing  process,  by  itself,  provides  enough  energy  to  generate 
sufficiently  strong  oscillations  post-failure  to  reduce  the  “effective  normal  stress/effective  friction”  to  the 
low  values  required  in  the  weak  core  region  near  the  cavity. 

We  are  implementing  an  approximate  representation  of  acoustic  fluidization  within  our  nonlinear  finite 
difference  code  and  have  applied  it  to  simulation  of  PILEDRIVER.  Two  different  schemes  are  being 
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implemented,  both  of  which  have  the  following  characteristics  in  common:  1)  Acoustic  fluidization  occurs 
only  after  failure  has  occurred  in  the  cell,  and  2)  Acoustic  fluidization  is  related  to  the  kinetic  energy  of 
each  cell.  So  the  basic  assumptions  of  the  calculation  are  that  once  a  cell  has  failed,  it  is  subject  to  acoustic 
fluidization,  and  that  the  acoustic  fluidization  occurs  only  after  the  kinetic  energy  reaches  some  critical 
threshold  level.  The  second  assumption  is  not  quite  right  -  acoustic  fluidization  is  related  to  acoustic 
vibrations  within  each  cell,  but  it  is  reasonable  to  expect  these  vibrations  to  be  approximately  related  to  the 
kinetic  energy  of  the  cell.  We  also  introduce  a  decay  constant  to  simulate  the  decaying  acoustic  vibrations 

d6  pp2 

in  each  cell  that  satisfies  the  equation  r  —  +  0  = - 

dt  2 

where  T  is  the  decay  time,  p  is  density,  v  is  particle  velocity,  and  9  is  the  kinetic  energy.  Acoustic 
fluidization  occurs  if  0  exceeds  a  threshold  level  in  the  rubbleized  cell.  The  effect  of  this  equation  has  been 
shown  in  preliminary  calculations  to:  1)  by  delaying  the  onset  of  9,  smooth  the  large  strength  reduction 
introduced  by  acoustic  fluidization  (and  to  smooth  the  velocity  pulses),  and  2)  slow  the  rapid  decrease  in  9 
at  the  end  of  the  step  function.  An  additional  modification  has  been  introduced  recently  to  allow  the  onset 
and  decay  times  to  occur  at  different  rates,  since  it  was  found  that  the  equation  above  causes  too  long  a 
delay  in  onset  time.  In  the  example  below,  we  discuss  only  the  case  for  x=0  which  corresponds  to  an 
instantaneous  onset  time  and  abrupt  stoppage  of  acoustic  fluidization.  Additional  calculations  with  varying 
onset  and  decay  rates  are  currently  being  conducted. 

Implementation  1  -  Failure  and  acoustic  fluidization  defined  by  damage  model 

First-principles  implementation  of  the  strong  vibrational  mode  of  the  acoustic  fluidization  mechanism  is 
beyond  the  capabilities  of  the  finite  difference  continuum  code  that  computes  the  ground  motion. 

However,  we  have  included  some  of  the  observed  phenomenological  features  of  the  acoustic  fluidization 
mechanism  in  our  computational  model.  Application  is  restricted  to  rock  elements  that  have  both 
previously  been  reduced  to  rubble  (using  the  Sammis  brittle  failure  model)  and  been  driven  to  a  sufficiently 
large  kinetic  energy.  These  rock  elements  are  then  assigned  lower  shear  strengths  through  the  use  of  a 
reduced  effective  pressure  or  an  increased  “acoustic  fluidization  pressure”  in  the  standard  friction  law. 

2 

We  implement  this  mechanism  in  terms  of  the  specific  kinetic  energy  density  (pv  /  2)  level,  called  Ek. 
When  this  kinetic  energy  goes  above  some  threshold  level,  Afonke,  acoustic  fluidization  occurs  and  the 
effective  pressure  consequently  drops  down.  In  the  finite  difference  code,  the  effective  pressure  field  is 
given  by  Pacf  =P-Pf ,  where  P  is  pressure  and  Pf  is  the  acoustic  fluidization  induced  pressure. 

The  stress  difference,  Y=b*  Pacf,  is  related  to  Pacf  using  a  friction  law  where  b  is  the  friction  coefficient. 

We  assume  that  acoustic  fluidization  induced  pressure  is  proportional  to  both  pressure  and  the  kinetic 
energy  density  in  the  mass  element  through  the  relation,  Pf=Ek/Afkemx  *P,  where  Afkemx  is  a 
normalization  factor.  The  effective  pressure  Pacf  would  drop  to  zero  if  the  kinetic  energy  is  greater  than 
afkemx.  Note  that  Pf  decreases  as  the  particle  velocity  and  kinetic  energy  in  a  cell  become  smaller  behind 
the  shock  front.  Thus,  the  strength  reduction  due  to  acoustic  fluidization  decreases  with  time  in  the  cell. 

This  section  contains  finite  difference  simulations  of  the  PILEDRIVER  explosions  that  were  made  to 
evaluate  both  the  Sammis  micro-mechanical  damage  mechanism  model  and  the  acoustic  fluidization 
mechanism  by  comparing  synthetics  and  observations.  We  start  first  with  a  series  of  runs,  aiming  to  find 
the  best  parameters,  afonke,  afkemx,  according  to  the  above  mechanism,  to  approximate  the  observations. 

Table  3  shows  results  such  as  failure  extent,  damage  extent,  acoustic  fluidization  extent,  final  cavity  size 
(Rc)  and  maximum  cavity  size  and  RDP  different  combinations  of  the  two  parameters,  afonke  and  afkemx. 
We  notice  that  the  final  scenarios  are  quite  sensitive  to  these  parameters,  and  generally  the  greater  value  of 
afonke  and  afkemx  corresponds  to  the  smaller  damage  zone  and  smaller  failure  zone.  We  find  that  the 
combination  Afonke=afkemx=le6  (erg/cm3)  (Case  Pa04  in  Table  3)  better  fits  the  cavity  size  and  RDP  for 
PILEDRIVER.  The  particle  motion  observations  and  synthetics  with  this  combination  of  parameters  at  four 
locations  are  shown  in  Figure  7. 
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Table  3  Parametric  studies  of  acoustic  fluidization. 


RUN-ID 

Afonke 

(erg/cm3) 

Afkemx 

(erg/cm3) 

Failure  extent 
(m)  (ifdam=4) 

Damage  extent 
(m)  (ifdam=l) 

Extent  of 
A.F.  (m) 

Rc  (final) 
(m) 

Rcmax 

(m) 

RDP.  phi 

(103  m3) 

Pa03 

le6 

le7 

367.6 

552.4 

212.2 

47.391 

48.31 

23 

Pa04 

le6 

le6 

436.8 

630.1 

212.2 

42.748 

56.95 

16 

Pa05 

le6 

5e6 

394.1 

590.1 

212.2 

50.045 

51.19 

26 

Pa06 

le6 

2e6 

436.8 

609.8 

220.6 

49.202 

54.98 

25 

Pa09 

le7 

le7 

367.6 

552.4 

129.0 

46.457 

47.19 

21 

PalO 

le5 

le7 

367.6 

552.4 

342.7 

47.480 

48.40 

22 

Pall 

2e5 

le6 

451.9 

630.1 

319.2 

43.090 

58.49 

17 

Figure  7.  Comparison  between  observations  (dashed  lines)  and  synthetics  (solid  lines)  made  with  acoustic 
fluidization.  afonke=afkemx=106  erg/cm3. 

The  simulated  waveforms  shown  in  Figure  7  match  some  of  the  characteristics  of  the  measured  waveforms 
but  are  not  as  good  a  fit  to  the  data  as  were  those  shown  earlier  in  Figure  1  (for  the  effective  stress  model) 
or  Figure  6  (for  the  Sammis  model  augmented  by  a  post-failure  shear  damage  model).  In  particular,  the 
very  simple  model  for  acoustic  fluidization  outlined  above  results  in  larger  negative  pulse  amplitudes  and 
shorter  negative  pulse  durations  than  those  measured.  Part  of  this  is  due  to  a  problem  in  the 
implementation  -  as  velocity  changes  from  positive  to  negative  the  kinetic  energy  reaches  zero  and  acoustic 
fluidization  abruptly  turns  off,  which  is  clearly  incorrect.  Introduction  of  a  decay  time  into  the  above 
acoustic  fluidization  algorithm  to  correct  this  problem  has  shown  the  potential  for  increasing  negative  pulse 
durations  and  bringing  the  acoustic  fluidization  results  into  better  agreement  with  the  measurements. 


However,  the  noise  in  the  waveforms  produced  using  the  Sammis  brittle  failure  model  (Figures  6  and  7)  is 
a  consequence  of  the  large  strength  reduction  for  failed  rock  compared  to  the  elastic  (infinite  strength) 
deviatoric  behavior  required  for  compatibility  with  the  model  in  the  damaged  rock  outside  of  it.  (See 
Rimer,  et  al,  1999,  for  a  discussion  of  the  effects  of  numerical  smoothing  of  the  discontinuity  between 
failed  and  damaged  rock  elements  for  the  PD  10  model.) 


Implementation  2  -  Failure  and  start  of  acoustic  fluidization  defined  by  laboratory  strength  model,  followed 
by  plastic  flow 


The  simulations  described  in  the  previous  section  are  noisy,  in  large  part  due  to  incompatibilities  both 
between  the  brittle  failure  model  and  possible  plastic  behavior  near  the  explosive  source  and  between  the 
elastic  pressure-volume  relation  for  undamaged  rock  and  the  nonlinear  loading  relation  measured  on  small 
samples  of  granite  in  the  laboratory.  In  order  to  incorporate  these  nonlinear  rock  behaviors,  we  removed  the 
Sammis  model  and  instead  used  the  laboratory  failure  surface  for  damaged  granite  as  the  trigger  point  for 
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acoustic  fluidization.  Calculations  were  run  with  this  failure  surface  and  the  polynomial  EOS  for  granite 
developed  from  compressive  loading  tests  on  Piledriver  rock  samples.  We  assume  plastic  yielding  and  a 
non-associated  flow  rule  (radial  return)  and  use  the  plastic  work  (epw)  done  on  a  rock  element  as  the 
threshold  to  determine  if  the  cell  is  rubbleized  (to  initiate  the  acoustic  fluidization  model).  The  acoustic 
fluidization  threshold  is  set  to  the  same  value  as  implementation  #1.  These  calculations  are  in  progress. 

CONCLUSIONS  AND  RECOMMENDATIONS 


IDG  is  in  the  process  of  digitizing  data  from  25  nuclear  tests  at  the  Degelen  test  site.  The  complete  data  set 
consists  of  81  near  field  waveforms  recorded  underground  at  shot  depth  and  124  (unique)  near  regional 
seismic  records.  IDG  has  digitized  27  near  field  waveforms  and  122  seismic  waveforms  (including  multiple 
records  of  the  same  signal),  and  will  complete  digitization  of  the  remainder  of  the  data  set  over  the  duration 
of  this  project.  The  data  set  of  near  field  and  regional  waveforms  will  be  delivered  to  the  Center  for 
Monitoring  Research  upon  completion  of  this  project. 

We  are  in  the  process  of  defining  the  optimum  procedures  to  model  this  data.  The  goal  is  to  develop 
material  models  that  are  consistent  with  the  data  and  have  a  realistic  physical  basis.  Work  to  date  at  SAIC 
has  focused  on  implementation  and  testing  of  improved  numerical  modeling  procedures  and  simulation  of 
near  regional  data.  We  are  testing  acoustic  fluidization  as  a  physical  mechanism  for  strength  reduction  in 
nonlinear  explosion  simulations.  Preliminary  results  are  promising,  but  additional  work  is  required  to 
improve  the  realism  of  the  models.  This  initial  study  has  focused  on  modeling  of  the  PILEDRIVER  data  set 
with  new  material  models.  During  the  next  year,  we  will  continue  this  modeling  effort  and  also  focus  on 
modeling  the  Degelen  data  set.  Near  field  data  will  be  modeled  with  our  nonlinear  finite  difference  code  as 
described  in  this  report  with  improved  material  models.  Near  regional  data  will  be  modeled  using  wave 
number  integration. 
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